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1 INTRODUCTION 



The alignment between the spin and velocity orientations of 
pulsars has been motivated both theoretically and observa- 
tionally. A number of mechanisms have been proposed that 
predict asymmetric kicks which lead to pulsar birth veloci- 
ties along the direction of their spin axis (e.g. Tademaru & 
Harisson 1975; Spruit & Phinney 1998; Cowsik 1998). More- 
over, pulsar-population syntheses based on those models 
predict certain correlations between pulsar observables and 
the degree of alignment that if present would lend support 
to the above mechanisms (e.g. Ng & Romani 2007; Wang , 
Lai k Han 2007; Kuranov, Popov & Postnov 2009). Obser- 
vationally, there is strong evidence for spin-velocity align- 
ment in individual pulsars, from X-ray imaging of pulsar- 
wind nebulae tori (e.g. the Crab and Vela pulsars; Caraveo 
& Mignani 1999; Ng & Romani 2004; Ng & Romani 2007), 
and in samples of pulsars where polarimetric observations 
can reveal the spin-axis orientation (Johnston et al. 2005; 
Rankin 2007; Noutsos et al. 2012). As was highlighted in 
Noutsos et al. (2012), without available information on the 
3D orientation of pulsar velocities and spin-axis directions, 



ABSTRACT 

This paper presents a detailed investigation of the dependence of pulsar spin-velocity 
alignment, which has been observed for a sample of 58 pulsars, on pulsar age. At first, 
our study considers only pulsar characteristic ages, resulting in no change in the degree 
of correlation as a function of age, up to at least 100 Myr. Subsequently, we consider 
a more reliable estimate of pulsar age, the kinematic age, assuming that pulsars are 
born near the Galactic plane. We derive kinematic ages for 52 pulsars, based on the 
measured pulsar proper motions and positions, by modelling the trajectory of the 
pulsars in a Galactic potential. The sample of 52 pulsar kinematic ages constitutes 
the largest number of independently estimated pulsar ages to date. Using only the 33 
most reliable kinematic ages from our simulations, we revisit the evolution of spin- 
velocity alignment, this time as a function of kinematic age. We find that the strong 
correlation seen in young pulsars is completely smeared out for pulsars with kinematic 
ages above 10 Myr, a length of time beyond which we expect the gravitational pull 
of the Galaxy to have a significant effect on the directions of pulsar velocities. In the 
discussion, we investigate the impact of large distance uncertainties on the reliability 
of the calculated kinematic ages. Furthermore, we present a detailed investigation 
of the implications of our revised pulsar ages for the braking-index and birth-period 
distributions. Finally, we discuss the predictions of various SN-kick mechanisms and 
their compatibility with our results. 
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statistical studies of pulsar spin-velocity alignment that are 
based on pulsar proper motions and polarimetric data have 
to rely on the projected spin-axis and velocity vectors onto 
the plane of the sky. Apart from the few exceptions men- 
tioned above, this is the case for all pulsars considered in 
Noutsos et al. (2012) and in the follow-up study presented 
here. 



Recently, Kuranov, Popov & Postnov (2009) used 
Monte Carlo (MC) population synthesis to examine the ef- 
fect of binary break-up and isolated pulsar kicks on the 
observed distribution of space-velocity magnitudes and the 
spin-velocity offset angles. The above authors compared the 
simulated distributions with those from the observed data 
sample of Rankin (2007) and Ng & Romani (2007), which 
helped them determine valid ranges for physical parameters, 
like the range of spin-kick-velocity offsets of a pulsar pop- 
ulation with various fractions of binary progenitors. How- 
ever, so far, there has been little attention to the depen- 
dence of the observed alignment on pulsar age. Indeed, if 
there is a mechanism that favours spin-velocity alignment, 
then one should expect to see such a correlation for young 
pulsars (< 1 Myr). For older pulsars (> 10 Myr), the Galac- 
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Figure 1. Probability distributions of the difference between the position angles of the pulsar spin (PAo) and velocity (PA V ), \£, for 
a range of characteristic-age intervals, r c ( m i n ) ^ r c < T c ( max ). Each distribution shown was generated from a large number of MC 
realisations of PAo and PA V , based on their published la errors. The number of spin-velocity angles considered in each distribution is 
shown as N psT , where the total number of pulsars (iVp Sr = 58) was considered for the distribution covering the entire investigated range 
of characteristic ages, i.e. ^ r c < 100, 000 kyr. The error bar shown at each 5°-bin of \I/ corresponds to the standard deviation of the 
bin's height across all MC iterations. The probability of rejecting the uniform ^ distribution, given the observed distribution, calculated 
with the Kolmogorov-Simirnov (KS) test, is shown as pks- The la confidence interval of pks, for each distribution of ^, was derived 
from the distribution of the KS statistic across the MC simulation. 



tic potential will alter the pulsar proper motions, diluting 
the effect and causing the correlation to weaken. The typ- 
ical timescale in which this happens can be estimated as 
[GMmw /(7rRM W 6h)]~ 1/2 - 10 Myr, where M MW 10 12 
M and Rmw = 25 kpc is the mass and radius of the Milky 



Way, respectively, and 5h = 200 pc is the Galactic plane 
thickness. Hereafter in this paper, we refer to the above es- 
timate as the "dynamical time" , tdyn • The effect of the galac- 
tic potential was perhaps in play in the studies of Johnston 
et al. (2005, 2007): the earlier study showed a clear spin- 
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velocity correlation, whereas the follow up work, using a 
separate, older sample of roughly equal size, showed no sig- 
nificant correlation. The critical difference between the two 
samples was the pulsar ages and, what is of equal impor- 
tance, the larger distances of the second sample compared 
to the first. For the majority of pulsars, their distance is 
estimated from their measured dispersion measures (DMs) 
and a model of the free-electron-density distribution. Us- 
ing the NE2001 electron-density model of Cordes & Lazio 
(2002) typically leads to uncertainties of the order of 20% on 
the distance (Cordes & Lazio 2003). Recent investigations 
by Gaensler et al. (2008) and Schnitzeler (2012) have shown 
that in a number of cases, e.g. for high-latitude pulsars, the 
discrepancy between the distance derived by the model and 
that obtained through parallax measurements can be even 
larger. Hence, corrections of pulsar proper motions for the 
Galactic differential rotation are typically less accurate for 
distant pulsars, which could, again, weaken an intrinsic cor- 
relation. 

The following sections of this paper investigate the ef- 
fect of pulsar age on the spin-velocity correlation, using the 
sample of Noutsos et al. (2012). Our sample is composed 
of 58 pulsars for which polarimetric and proper-motion in- 
formation was available; in a few cases (i.e. J0534+2200, 
J1709-4429 and J1952+3252) the spin-axis direction was 
derived from fitting for the orientation of the pulsar wind 
nebula tori in X-rays. Detailed information about the ob- 
servations, the data selection and reduction can be found in 
Noutsos et al. (2012) and the references therein. The statis- 
tical analysis performed on the total sample in our previous 
paper showed that spin-velocity alignment is favoured at the 
99% confidence level. Furthermore, using toy-model simula- 
tions, we showed that the distribution of spin-velocity align- 
ment angles (considering polarisation and orthogonal mode 
ambiguities) resembles more a distribution of truly aligned 
configurations and less that of those where spin and velocity 
are orthogonal. Finally, it was shown that the conclusions 
are not biased by systematics, e.g. significantly aligned sub- 
sets of data or the choice of PAo . 

In this paper, we examine the dependence of spin- 
velocity alignment on pulsar age. As a first step, we take 
the spin-down ages of our pulsar sample (r c = P/(2P)) at 
face value and examine the distributions of the offset angle 
between spin and velocity, ^, as a function of r c . As a fol- 
lowing step, we attempt to provide an alternative measure 
of the pulsar ages by assuming that our pulsars were born 
close to the Galactic plane (GP): we trace their Galactic 
trajectories back in time, using their current proper motions 
and positions and calculate their kinematic ages, tkin, as the 
time interval between their current position and past inter- 
sections with the GP. For the sample of pulsars, for which 
the determination of tkin was possible, we revisit the de- 
pendence of spin-velocity alignment on age using the more 
reliable tkin instead of r c . Finally, we investigate what effect 
large uncertainties on pulsar distance have on the determi- 
nation of kinematic ages and the spin- velocity correlation. 
Also, the assumptions that pulsars are born with spin peri- 
ods, Po, that are much shorter than those observed today, 
and that they spin down by converting their rotational en- 
ergy to pure magnetic-dipole emission, i.e. P oc 1/P, is also 
tested for our sample of pulsars. 




100,000 



Figure 2. Scatter plot of the probability of rejecting the uniform 
distribution as a function of characteristic age (r c ), for the \& 
distributions along the main diagonal in Fig. [T] The data points 
shown correspond to the values of pks f° r the distributions along 
the main diagonal in Fig. ^ The number of pulsars considered 
in the corresponding \I/ distribution is also shown at the position 
of each marker. The length of time corresponding roughly to the 
dynamical time of ~ 10 Myr is shown with a thick, vertical grey 
line (see Section [TJ. 



2 DATA ANALYSIS 

The statistical procedure followed in this paper is similar 
to the one of Noutsos et al. (2012). Firstly, based on the 
measurement uncertainties on the position angle of the sky- 
projected pulsar velocity, PA^, and that of the sky-projected 
spin axis, PAo, we generated a large number of MC data sets, 
each having equal size to the original sample of 58 pulsars. 
The MC data sets were generated by randomly drawing val- 
ues from 58 Gaussian distributions with means and standard 
deviations equal to the published PA values and their uncer- 
tainties, respectively. Secondly, for each of the data sets we 
calculated the values of ^ = PAo — PA V , while considering 
the polarisation ambiguity (PAo =b 180°) and the ambiguity 
due to orthogonal- mode emission (PAo =b 90°). Due to these 
ambiguities, the minimum difference between PAo and PA V 
can be defined inside the interval, —45° ^ ^ ^ 45°; but 
since we are only interested in the absolute offset between 
the spin and velocity vectors, we only considered the ab- 
solute value of the above difference. Hence, all values of ^ 
were distributed inside the interval, 0° ^ ^ ^ 45°. During 
the last step, we quantified the degree of spin-velocity cor- 
relation, given the observed ^ distributions, using the one- 
sample Kolmogorov-Smirnov (KS) test. This test provides, 
via its test statistic, D, the probability that our observed 
data set of ^ values is drawn from an arbitrary theoreti- 
cal distribution, e.g. the uniform distribution, p($?) = 1/45 
deg -1 . Throughout this paper, we quote the probability of 
rejecting the uniform distribution under the KS test, pks' 
this is defined as 1 minus the probability that the observed 
^ distribution is drawn from a uniform one. In the final step, 
we generated the cumulative histograms of all MC data sets 
by summing the values in each bin and calculated the value 
of pks corresponding to the mean value of D from all MC 
runs. 

We examined the distributions of ^ and the correspond- 
ing pks values for different subsets of our pulsar sample, 
which were defined according to intervals in r c . The consid- 
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ered range of characteristic ages was r c < 100 Myr, which 
contained the entire pulsar sample. The range was split into 
intervals with sizes ranging from one to four orders of mag- 
nitude in t c . Fig. [1] shows the tabulated distributions of ^, 
each distribution corresponding to that subset of pulsars 
which have r c in the interval whose lower and upper bounds 
are shown in kyr along the first column and top row, respec- 
tively. The number of pulsars belonging to each age interval 
is also quoted with each distribution. In addition, to visu- 
alise a possible trend in pks with increasing r c , we show in 
Fig. [2] the pks values of the distributions along the main 
diagonal of Fig. [1] 

At first glance, all ^ distributions appear non-uniform 
with most of the ^ values distributed below 20°. In general, 
the incremental addition of older pulsars is accompanied by 
increasing pks values, i.e. pks(t c < 100 kyr) < pks(t c < 
1,000 kyr) < p K s(r c < 10,000 kyr), etc. The distribution 
that was generated from all pulsars in our sample shows the 
highest probability of rejecting uniformity, with the KS test 
reporting pks ~ 99.8%. However, the probabilities corre- 
sponding to discrete intervals of r c , as is shown in Fig. [2j 
show no correlation between £>ks and r c . Interestingly, we 
also have £>ks(t c ^ kyr) > £>ks(t c ^ 100 kyr) > £>ks(t c ^ 
1,000 kyr), etc. Overall, there is no evident trend towards 
more uniform distributions of ^ (i.e. smaller values of £>ks) 
with increasing characteristic age. Nevertheless, it is true 
that the histograms for the very young pulsars (< 100 kyr), 
those with ages between 100 and 1,000 kyr, as well as the 
very old pulsars 10, 000 kyr) are derived from a small 
sample of ~ 10 pulsars. Therefore, our confidence in the re- 
ported probabilities for those cases is generally lower than 
for the rest of the distributions, which is also reflected in the 
size of the errors on pks: these were derived from the 68% 
confidence limits (CLs) of the KS statistic, D. The excep- 
tion to this is the distribution for 1, 000 ^ r c < 10, 000 kyr, 
which is based on a sizeable sample of 32 pulsars but has 
also a large lower uncertainty on pks (i-e. 23% lower error). 

In conclusion, if we take the above results at face value, 
we can exclude (at the 90% confidence level) that spin and 
velocity orientations are completely uncorrelated, for the en- 
tire investigated range of characteristic ages. Yet, if there is 
an age-dependent smearing of the spin-velocity alignment, 
it may only become effective beyond ~ 10 8 y, assuming r c 
to be reliable (but read Section [3}. However, as was noted in 
the introduction, this is hard to reconcile with the dynami- 
cal time of ~ 10 Myr, being an order of magnitude shorter 
that this limit. 



3 PULSAR KINEMATIC AGES 

Pulsar age is a critical parameter in our investigation of 
whether pulsar spin vectors tend to be aligned with their 
velocities for young pulsars. However, except for the few 
cases where historical pulsar-supernovae associations can 
provide good estimates of pulsar ages, the only other mea- 
sure we possess of a pulsar's age is its spin-down or char- 
acteristic age (t c oc P/P). Unfortunately, in many cases it 
can be readily seen that the characteristic age is an unreli- 
able estimate of the true age, ttme'- e.g. tracing the outward 
proper motions of PSR B1951+32 and PSR B0538+2817 
from the centres of their respective supernova remnants 



(SNRs) has revealed that their true ages are £ t rue ~ 0.5r c 
and ~ 0.05t c , respectively (Migliazzo et al. 2002; Kramer et 
al. 2003). This is largely due to the implicit assumptions of 
the above expression for r c , which are that pulsars are born 
with periods much shorter than those observed today and 
that they spin down by converting their rotational energy to 
low-frequency magnetic-dipole radiation, while maintaining 
a constant magnetic- field. An additional source of system- 
atic uncertainty in the determination of r c comes from the 
"Shklovskii effect" , whereby the apparent spin-down rate, P, 
is contaminated with kinematic contributions like the rela- 
tive transverse velocity between the pulsar and the solar sys- 
tem barycentre (Camilo et al. 1994); the Galactic differential 
acceleration and vertical acceleration introduce additional P 
components (Damour & Taylor 1991). The Shklovskii effect 
on the measured P, and hence on r c , becomes important for 
high-velocity, nearby pulsars, whereas Galactic differential 
acceleration and vertical acceleration are important for dis- 
tant and high latitude pulsars, respectively. In general, the 
contribution of these dynamical effects to the value of P is 
negligible for slow-spinning, non-recycled pulsars. However, 
for nearby, high- velocity millisecond pulsars, with values of 
P that are ~ 5 orders of magnitude smaller, the above ef- 
fects may account for a significant fraction of the measured 
P (see e.g. Section 8.2.4 of Lorimer & Kramer 2005). 

Given the aforementioned uncertainties associated with 
the determination of pulsar ages based on spin down, it 
is important to explore alternative methods of determin- 
ing those ages that are independent of a spin-down model. 
To date, there are only a few non-recycled pulsars whose 
ages have been determined through independent methods. 
Some of them exhibit a large discrepancy, up to an order 
of magnitude, between r c and the independent estimates. 
The discrepancy can reveal r c as being either an over- or 
under-estimate of the true age of a pulsar. For instance, 
the true age of PSR J0538+2817 is roughly an order of 
magnitude younger that its characteristic age (Kramer et 
al. 2003). Conversely, the westward proper motion of PSR 
J 180 1—2451 suggests that the pulsar's true age must be sig- 
nificantly higher than its characteristic age, if the pulsar 
is associated with in the SNR G5.4— 1.2 (Gaensler & Frail 
2000; also see Zeiger et al. 2008). Furthermore, the outwards 
proper motion of PSR J 1932+3252 from the geometric cen- 
tre of its parent SNR CTB 80 implies a true age of ~ 51 
kyr for this pulsar (Zeiger et al. 2008). However, these ex- 
amples are of typically young pulsars with known supernova 
associations; it is quite possible that a large number of older 
pulsars, for which we have no independent clues about their 
ages, have significantly different true ages to their charac- 
teristic ones. 

Under the assumption that pulsars spin down entirely 
due to electromagnetic radiation, it is generally assumed 
that the spin-down rate is proportional to some power law 
of the spin rate: this is expressed as 

v = -Kv n , (1) 

where v — 1/P is the spin frequency; K is a proportionality 
factor, often assumed to be constant, that depends on the 
physical properties of the pulsar, i.e. the moment of iner- 
tia, the magnetic-dipole moment and the inclination angle 
between the spin- and magnetic axes; n is the pulsar brak- 
ing index (Manchester & Taylor 1977; Shapiro & Teukolsky 
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1983). The general solution of the above differential equation 
for the true age of a pulsar, ttrue, is 



(n-l)P 



if n ^ 1 (2) 



Here, we have implicitly assumed that apart from K, the 
braking index also remains constant with time. However, 
the observed trends in a number of pulsar properties with 
age, e.g. the magnetic inclination angles, suggest that K 
and n are functions of time, implying that pulsar magnetic 
fields decay with time and/or that their magnetic axes tend 
towards alignment with their spins, as they age (Tauris h, 
Manchester 1998; Young et al. 2010). If this is true, then the 
above solution is invalid; in that case, Eq. [1] is generalised 
as 



m = -K(t)Ht)] n(t \ 



(3) 



which if solved would give the true age of a pulsar but whose 
form is unknown and likely varies between pulsars. Accord- 
ing to Vranesevic & Melrose (2011), if the evolutionary mo- 
tion of pulsars on the P-P diagram constitutes a current 
that is a conserved quantity, i.e. the observed distribution is 
maintained in time, then pulsar spin-down can be generally 
expressed as P = —dC(P)/dP, where C(P) is the unknown 
potential responsible for the pulsar current. In the specific 
case of pulsar evolution under a constant n, the above gen- 
eralised spin-down law reduces to Eq. [TJ 

Under the usual assumptions of Po/P = and n = 3, 
Eq. [2] simply reduces to the characteristic age or spin-down 
age of a pulsar 

P 

Tc = • (4) 

2P 

However, all reliably measured braking indices to date are 
less than 3; amongst them there are extreme cases where 
n « 1, like that of the Vela pulsar and PSR J1734-3333 
(Espinoza et al. 2011). Hence, if these pulsars were born 
fast and have evolved under a constant braking index, then 
t c would be an underestimate of their true age. The age of 
pulsars with n = 2 would then be closer to 



T 2 - 2t c 



(5) 



where we have assumed again that Po/P = in Eq. [2] Fur- 
thermore, for pulsars that have evolved with n — 1 even T2 
would be an underestimate. For this special case, Eq. [T]has 
a special solution, i.e. 



tt 



2r c In ( ^~ 

Po 



if n ■ 



(6) 



which cannot be simplified further. Nevertheless, we can ar- 
bitrarily choose Po to be very small as to have again Po/P ^ 
for non-recycled pulsars. Choosing Po = 1 ms, close to 
the break-up limit for a neutron star (e.g. Chakrabarty et 
al. 2003), we then have 



r\ — 2t c In 



1 ms 



(7) 



where P is in ms. 

Finally, it is easy to show that pulsar age, as defined 
by Eq. [2] increases monotonically with decreasing n, for any 
value of < Po < P and any n (this also remains true 
around the special case of n — 1). Consequently, for every 



pulsar, n > T2 > t c : i.e. larger values of n mean that the 
pulsar spins down more rapidly. 



3.1 Characteristic Age Bias due to P and n 

Invoking Eqs. 11121 and HI the ratio between ttrue and r c , at 
an arbitrary time, t, can be expressed as a function of the 
ratio, n, measured at the present time, t±: 



s,\ ttrue 



n 



i + x/ti 
i + x/t 



(8) 



where x = Po _1 /[K( n ~ !)]• Since Po is usually unknown, 
we can express the constant \ as a function of r\ and t\ 
by solving Eq. ^ and replacing the expression for P(t) in 
t c = P/(2P): 



"(*) 



1 + (p/ri - 1) 



(9) 



where we have set p = 2/(n — 1). 

As an example, we can use Eq. [9] to predict the ex- 
pected discrepancy between ttrue and r c , for a particular pul- 
sar, in future. As was mentioned earlier, PSR J0538+2817 
has t c = 600 kyr and ti = 30 kyr, which translates to 
r\ — 0.05. Assuming n — 3, according to Eq. [9] the ratio 
ttrue/ V c becomes r ~ 0.95, when ttrue = 10 Myr. This result 
is expected, because, even though today's period of PSR 
J0538+2817 is comparable to its birth period, as the pulsar 
becomes older the assumption that Po <C P becomes more 
valid. Moreover, if n = 2, Eq. [9] gives r = 1.8, for t = 10 
Myr. In the more extreme case of the Vela pulsar, the brak- 
ing index could be as low as 1.2 (la lower bound), which 
would mean r = 6.3att=10 Myr. Finally, for the special 
case where n — 1, we simply have r(t) = ri(t/ti); for PSR 
J0538+2817, it would mean that r = 16.7. So, the age bias 
in t c due to the unknown birth period of pulsars is less im- 
portant for older pulsars but can still be as big as an order 
of magnitude for young pulsars (e.g. PSR J0538+2817). On 
the other hand, the bias due to the braking index is far more 
important for old pulsars. If several of the old pulsars in our 
sample (> 10, 000 kyr) have evolved under a constant brak- 
ing index that has been close to 1, like that of the Vela pulsar 
at present, then their characteristic ages could be as much 
as 10 times younger than their true ages (see Section [472]) . 

In the worst case, both of the above effects could be in 
operation for the pulsars in our sample, which would mean 
that if we replaced r c with ttrue, a number of pulsars would 
belong to different age intervals in Fig. [T] More specifically, 
middle-aged pulsars according to r c (=100 — 10,000 kyr) 
would be shifted to the distributions for ttrue < 100 kyr, 
due to their large Po. Similarly, some of those middle-aged 
pulsars would belong to the distribution for ttrue > 10, 000 
kyr, due to n < 3. In the following paragraphs we consider 
an alternative method of estimating the age of the pulsars 
in our sample and re-examine the spin-velocity alignment 
as a function of age, based on those new estimates. 

3.2 Methodology 

It is generally believed that short-lived, massive OB stars, 
that dwell near the GP, are the likely progenitors of pul- 
sars (e.g. Amnuel, Guseinov & Kustamov 1986; but see also 
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Blauw 1985). It is therefore a reasonable assumption that 
pulsars are born at Galactic heights of |zbirth| <C Rmw- 
Based on this assumption, pulsar ages can be estimated us- 
ing a kinematic analysis, by calculating the length of time 
required for a pulsar to travel between its location at birth 
and that at present, through the gravitational potential of 
the Galaxy. Such a method has been discussed in the litera- 
ture by various investigators (e.g. Lyne, Anderson & Salter 
1982; Brisken et al. 2003). In general, pulsar ages derived 
from such kinematic analyses are free of the biases discussed 
in the previous section and, thus, can be more reliable than 
characteristic ages. In the following paragraphs, we describe 
the method used in this paper to derive "kinematic ages", 
hereafter tkin, for the pulsars in our sample. Our goal is 
to re-examine the dependence of spin-velocity alignment on 
pulsar age, this time using the more reliable kinematic ages. 

3.2.1 The pulsar trajectories through the Galaxy 

The calculation of a pulsar's path through the Galaxy, 
e.g. between its birth site and its present position, requires 
a numerical integration of the equations of motion through 
a Galactic gravitational potential. Primarily, this calcula- 
tion requires knowledge of the following pulsar observables, 
which form the initial conditions for the integration, at 
present time: 

- The pulsar's position in Galactic coordinates, (/,£>). 

- The pulsar's distance, d. 

- The pulsar's proper motion, (/ij,/ib). 

Proper motions have been measured for hundreds of pulsars 
and can be found in the literature (e.g. Hobbs et al. 2005). 
On the other hand, pulsar distances are only accurately 
known for a few tens of pulsars, from parallax measurements 
and SNR associations. Nevertheless, the distances of most 
pulsars can be crudely estimated via their dispersion mea- 
sures (DM) and a model of the Galactic free-electron density 
distribution (Cordes & Lazio 2002; Schnitzeler 2012). 

In addition, the numerical integration requires an as- 
sumption about the Galactic gravitational potential. Several 
potentials exist in the literature (e.g. Kenyon et al. 2008; 
Kuijken & Gilmore 1989). In this work, we have chosen the 
potential of Paczynski (1990). These potentials are generally 
expressed in right-handed Galactocentric Cartesian coordi- 
nates (GCC), (x, y, z), where the x axis is directed from the 
Sun to the Galactic centre (GC). Since the pulsar parame- 
ters are always measured with respect to the solar system 
barycentre (SSB), it is convenient to express them also in 
GCC. The pulsar's present velocity in GCC, as measured 
with respect to the solar system barycentre (SSB), is then 
given by 

v x =v r sin / cos b + fiid cos I — /j.bdsml sin b 

v y — — v r cos / cos b + \iid sin / + \±\>d cos / sin b (10) 

v z —v r sin b + /ifedcos b 

The initial conditions for the numerical integration criti- 
cally depend on the unknown radial velocity, i.e. ussb^i-)- 
Like with Zbirth, we have also numerically explored a range 
of possible v r values. The choice of this range was moti- 
vated by the work of Hobbs et al. (2005), who compiled a 
large catalogue of 233 pulsar proper motions, derived from 



timing-parallax measurements. Their work showed that the 
distribution of observed magnitudes of pulsar transverse ve- 
locities, along one dimension (ID distribution), is roughly 
contained within 1,000 km s _1 . Hence, in this work we have 
explored the range \v r \ ^ 500 km s _1 , where positive v r is 
defined as being directed away from the observer. As before, 
we defined a fine grid with Av r = 1 km s _1 and calculated 
the pulsar's trajectory for each of the grid values. 

We can now numerically perform the integration 
of the equations of motion from the present values of 
TGCc(x,y, z) and ussbOe, 2/, 2) to the values at birth, 

r , GCc(^birth, 2/birth, ^birth) and ^SSB^birth, 2/birth, ^birth)- 

As was mentioned earlier, we assume that pulsars are 
born very close to the Galactic plane. In this work, we 
have numerically explored a range of birth heights, with 
l^birthl ^ 100 pc. The chosen range was motivated by the 
work of Reed (2000), who used the observed latitudes and 
magnitudes of OB stars to determine the scale height of 
their distribution to /iob ~ 45 pc. In his work, Reed as- 
sumed that the optical extinction of the Galactic disc in the 
B-band, within 100 pc, is 2 magnitudes kpc -1 . In order to 
numerically simulate a continuous sampling of the chosen 
range in Zbirth, we formed a very fine grid of birth heights 
with A^birth = 0.2 pc. Each intersection of the pulsar's tra- 
jectory with grid elements was considered a possible birth 
site. The birth location along the other two coordinates, 
i.e. (xbirth, 2/birth), was left unconstrained. 

Finally, a number pulsars, e.g. PSR J1900-2600 (r c ~ 
100 Myr), may have crossed the Galactic plane several times 
in their lifetime before reaching their current position. If un- 
restricted, and depending on the value of v r , our simulation 
would allow pulsars that are bound to the Galactic gravita- 
tional potential to cross the GP indefinitely. Therefore, we 
decided to place a strict upper limit on the length of our 
integration. In particular, for each pulsar, we considered all 
possible intersections with the GP within r±. Hence, our sim- 
ulation places an artificial upper limit on the kinematic age 
of each pulsar, i.e. tkin ^ n. This limit corresponds to the 
case in which the pulsars are born with Po — 1 ms and spin 
down exponentially on a characteristic timescale of T2 = 2r c 
(see Eq. [7]). Given the spin periods of the 58 pulsars in our 
sample, this limit ranges between n ~ 7r c and 15r c , with 
52 out of the 58 pulsars having n > 10r c . Therefore, we 
deemed n a long-enough time scale for the purpose of our 
simulation. In Section 13.2.51 we discuss in more detail the 
implications of having multiple intersections with the GP. 

Following the formation of the above grids, for each pul- 
sar we calculated tkin for all (v r , ^birth) pairs that produced 
one or more intersections within t\. Up to this point, our 
method does not consider the actual probability distribu- 
tions of v r and Zbirth- In reality, those distributions are far 
from uniform, as it is more likely that a pulsar is born closer 
to the plane than away from it; similarly, it is more likely 
that a pulsar is travelling at ~ 100 km s _1 than at 1, 000 km 
s _1 . At the next step of our simulation, we consider those 
two probability distributions. 



3.2.2 The distribution of birth heights 

The work by Reed (2000) has shown that the probability 
density of the birth heights of OB stars, the likely progeni- 
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Figure 3. Posterior PDFs derived from our MC analysis for (a) PSR J0454+5543 and (b) PSR J0814+7429. (Central panels) The poste- 
rior probability, ^(vbirthj £kin)> from a large number of MC realisations of z and ^birth? randomly drawn from the exponential distribution 
of Reed (2000) and the Maxwellian distribution of Hobbs et al. (2005), respectively. (Bottom panels) The marginalised posterior prob- 
ability, p^birth)? derived from p(^birth^kin) by integration over all t^ n . (Right panels) The marginalised posterior probability, p(tkin) 5 
normalised by r c , derived from p(fbirth 5 ^kin) by integration over all ^birth- For PSR J0454+554, the simple form of p(tkin) allows us to 
quote the 68% (dotted line) and 95% (dashed-dotted line) confidence intervals around the peak value (solid line). For PSR J0814+7429, 
the p(tkin) nas a complex shape with multiple peaks of comparable probability. Hence, we cannot confidently quote a value for tkin5 m 
this case, only the age at PDF maximum is reported, also stating that this is an unreliable age estimate (see Table 



tors of pulsars, falls exponentially with z: i.e. 

m(z) ocexp(-H//i B) (11) 

where the exponential scale height, Hob was determined in 
that work to be 45 pc. In order to incorporate this prob- 
ability distribution into our simulation, we first generated 
a large number of values of z from 9l(z), allowing only 
\z\ < 100 pc. It is easy to show that the considered range 
of birth heights covers w 90% of d\(z). For each generated 
z value, we found its nearest neighbour in the grid of Zbirth 
and kept all the corresponding values of tkin that had been 
previously calculated for that Zbirth- After doing this, we are 
left with a large number of (v r , Zbirth) pairs and correspond- 
ing tkin values, appropriately weighted based on Eq. [11] 

3.2.3 The distribution of birth velocities 

In the previous paragraph, the consideration of the distribu- 
tion of birth heights ensured that smaller values of |zbirth| are 
weighted higher than larger ones. However, up to this point 
of our simulation all pulsar velocities, \v r \ ^ 500 km s _1 , are 
treated as equally probable. The distribution of pulsar ra- 
dial velocities, v r , is unknown, so it is not possible to directly 
follow the previous procedure, as was done for Zbirth- Never- 
theless, based on the proper motions of 233 pulsars, Hobbs et 
al. (2005) conclude that the distribution of pulsar velocities 
at birth, ^birth, can be fitted with a single Maxwellian: 

^) = \/|4exp(-^) (12) 

where a = 265 km s _1 is the standard deviation of the 
observed ID pulsar velocities, in either the longitudinal or 
latitudinal direction. Before we can use i^(vbirth) with our 
data, we need to transform the sample of v r values, from 
the previous step, to ^birth- This transformation entails the 
following two steps: 



- Correction of vssb for the motion of the local standard 
of rest (LSR), which affects the presently observed pulsar 
velocities. Hence, for each value of v r , we subtracted the 
velocity of the LSR, V©, which we assumed to be azimuthal 
at the Solar circle (Galactocentric distance of Rq = 8.5 kpc), 
directed clockwise, seen from the North: i.e. 

vlsr(x, y,z) = vssb(x,2/, z) + V®y (13) 

where we have chosen V© = —225 km s _1 , following the 
value used in Harrison et al. (1993). 

- Correction of the pulsar's velocity at the birth loca- 
tion, (xbirth, 2/birth, ^birth), for the in-situ Galactic rotation, 
V r Mw(^birth, 2/birth). This calculation requires first the deter- 
mination of the pulsar's birth location, which was performed 
by numerical integration of the equations of motion using 
vlsr- The resulting velocity is the pulsar's velocity at birth, 
relative to the local ISM: 

Vbirth = ^LSR^birth, 2/birth) + ^Mw(^birth, 2/birth) (14) 

whose magnitude, ^birth, is always positive. It should be 
noted that in the above calculations we have assumed a flat 
rotation curve for the Galaxy, i.e. no radial or vertical de- 
pendence of Vmw • 

3.2.4 The posterior probability distributions 

The above procedure transforms the original (v r , ^birth) pairs 
to pairs of (vbirth, 2birth)- In the last step, we first divided 
the ranges of v hirth e [vwh^wh] and t kin G [CnVSn*], 
for each pulsar, into 1,000 bins: on the ^birth-tkin plane, 
these form a 2D grid with bin size, £i>birth x £tkin, where 

c 1 / max min \ j c, 1 //.max <min\ 

O^birth — 1,000 V^birth — ^birth ) ana Otkin — 1 000 (Jkin — l kin )• 

Next, we generated a large number of values of v from $)(v). 
Note that the range of ^birth, although implicitly restricted 
by the chosen range of v r , can vary significantly between pul- 
sars with different locations and proper motions. For each 
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of the randomly drawn values of v, we calculated the num- 
ber density of the data from our simulation that fell into 
each of the 2D rectangular bins. The resulting 2D density 
plot constitutes the 2-dimensional (2D), posterior probabil- 
ity density function (PDF): 

/+100 pc 
dZbirth X 
-100 pc 

X £(^birth, tkin)^(^birth)^H(^birth), (15) 

where £(^birth, tkin) is the likelihood function that depends 
implicitly on the model parameters, v r and Zbirth — becom- 
ing zero outside their investigated ranges — and H(x) is the 
Heaviside function that simply indicates that all solutions 
corresponding to tkin > t\ are rejected. Fig. [3] shows two 
examples of grayscale density maps of p(^birth, tkin), where 
tkin is normalised by r c . Notice that we derived p(vbirth, tkin) 
by marginalising Zbirth over the considered range. 

A direct product of the above simulation were the 
marginalised posterior PDFs for ^birth and tkin, 

poo 

P(^birth) = / P(^birth, £kin)dtkin 
J 
/* oc 

P(tkin) = / P(^birth, t k in)^birth (16) 
J 

which are shown in the bottom and right-hand-side panels 
of Fig. [3] respectively. 

Evidently, despite the well-defined form of i^(vbirth) 
(i.e. the single- mode Maxwellian of Hobbs et al.), the shape 
of p(^birth) differs significantly from the prior distribution. 
The main reason for this is the multiplicity of solutions for 
tkin, for a given i>birth, which arises from the parametric pair 
of equations that relates ^birth and tki n : i.e. 

^birth = ^birth(^birth, V r ) 

tkin = tkin(^birth, V r ) (17) 

And, of course, the same is true for p(tkin), which means 
that we cannot recover the intrinsic PDF of tkin- It is im- 
portant to clarify that the multiple solutions to which we 
are referring, here, do not arise from multiple intersections 
with the GP but are the result of the multiple possibilities 
for (fbirth, ^birth) that lead to the same tkin, and, vice versa, 
the multiple (tkin, ^birth) that lead to the same ^birth- Al- 
lowing for multiple intersections adds an additional degree 
of complexity to the derived solutions, which is discussed in 
Section 1X231 

Finally, it should be stressed that, given the lack of 
prior information, e.g. association of the birth location with 
a spiral arm or high-density region, all valid solutions of 
Eq. [TTlin our simulation were equally weighted. 

3.2.5 Multiple Intersections with the Galactic Plane 

For the majority of pulsars in our sample, our simulation 
produces only a single intersection with the GP, within n 
and for all investigated v r values. In contrast, Fig.[4ji shows 
an example trajectory for PSR J1900— 2600 (black line), be- 
tween the pulsar's present position and that n ~ 608 Myr 
ago, if v r — —250 km s _1 . Under these particular assump- 
tions, it is seen that the pulsar crosses the GP 6 times, ex- 
iting the plane in between, which would thus produce 6 dis- 
crete solutions for tkin- The p(tkin) for this pulsar is shown 



in Fig. [4)3, where one can see that it is composed of mul- 
tiple peaks. However, it should be cautioned that the cor- 
respondence of the discrete peaks contained in p(tkin) to 
the crossings of the GP for this particular trajectory is not 
necessarily one-to-one, as the former PDF is the result of 
marginalisation over all v r values and, hence, over all pos- 
sible trajectories; these may produce a different number of 
intersections and at different times. Here, we clarify that our 
simulation treated each possible solution for each trajectory 
up to ti as equally likely. Nevertheless, the different prob- 
abilities of each peak in the PDFs arise from the following: 
(a) as was stated above, certain v r values may lead to fewer 
number of intersections than others, which when marginalis- 
ing over v r would lead to lower weighting at those tkin values 
corresponding to the least occurring intersection number. 
This can be made clearer in Fig. where the trajectory 
for v r — —250 km s _1 (black line) is compared to that for 
v r — 150 km s _1 (grey line). The latter produces an addi- 
tional intersection, i.e. 7, within n, compared to the former 
case, which contributes additional solutions towards p(tkin) 
at roughly the same tkin (note that tkin at intersections 6 
and 7 is roughly the same on a log scale — see below), (b) 
In addition, for a single v r , the crossing duration of an in- 
tersection (i.e. the time a pulsar spends in \z\ ^ 0.1 kpc), 
as well as the time interval between intersections affects the 
width and height of the peaks in p(tkin)- In particular, if 
e.g. two successive intersections occur close in time (i.e. a 
small interval between peaks) , then depending on the width 
of the peaks (i.e. the crossing duration) the corresponding 
probability at a certain tkin may be the sum of probabili- 
ties due to two or more intersections. And, of course, since 
p(tkin) is the average over all values of v r , the shape of the 
peaks in the PDF is also determined by the spread in tkin 
corresponding to each intersection, which is caused by the 
range of travel times for the different v r values, (c) Finally, 
the PDF of tkin is presented on a log scale with constant 
bin size, A log(tki n /Myr) = const. As a consequence, uni- 
formly distributed probability densities as a function of age 
on a linear scale would correspond to non-uniform probabil- 
ity densities on a logarithmic scale: this is simply expressed 
as Ap(t k in)/A log tkin oc tki n Ap(tkin)/Atkin- In other words, 
the number of solutions per age interval would statistically 
increase with tkin, favouring higher probabilities at older 
ages; however, this effect is mainly evident in PDFs having 
multiple peaks separated by at least an order of magnitude 
in tki n : as we explain in the following section, pulsars having 
such PDFs have been excluded from further analysis. 



3.3 Results 

3.3.1 The determined kinematic ages 

The majority of pulsars in our sample, i.e. 33 out of the 52, 
have a p(tkin) that exhibits a single, clear peak inside the 
investigated range (e.g. Fig.[3ji). Such cases allow us to con- 
fidently assign the most probable kinematic age to each pul- 
sar (i.e. the peak value) and quote the associated 68% upper 
and lower CLs around the peak. However, for the remaining 
19 cases, the corresponding p(tkin) is complex, having mul- 
tiple peaks — in many cases, of comparable significance. An 
extreme case of this category is the PDF of PSR J0814+7429 
(see Fig.[3]3), for which n ^ 1.7 Gyr, during which this pul- 
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Figure 4. (a) Trajectory of PSR J1900— 2600 between the pulsar's present position (circled dot) and that at t\ ~ 608 Myr in the 
past (cross), assuming v r = —250 km s _1 (black curve) and v r = 150 km s _1 (grey curve). The trajectories are projected on the y-z 
plane that is perpendicular to the sightline between the Sun and the GC. The intersections with the GP along the two trajectories are 
numbered in ascending chronological order, the most recent being the first (1). (b) PDF of tkin for PSR J1900— 2600, as was derived from 
our kinematics simulation for tkin ^ r\ (grey-filled histogram); for completeness, the histogram of solutions beyond r\ (marked with a 
dashed line) is shown as a grey outline up to 1 Gyr. 




y(kpc) "bi rth ( km s " 1 ) 



Figure 5. (a) Grey curve: 40-Myr trajectory of PSR J1932+1059 for the case of v r = —250 km s -1 , shown on the y-z plane (see 
Section 13.21 for a definition). The pulsar's current position is shown with an open, black circle, where the yz component of its present 
velocity assuming the above v r is also shown. The pulsar positions along the path at r c , T2 and n ago are also indicated with black 
bullets. Pulsar positions beyond n, here shown with a grey, dashed line, were not considered in our study. The horizontal, black line 
[z = kpc) and surrounding grey lines (0 < \z\ ^ 0.1 kpc) mark the range of possible ^birth values considered in our study; the range of 
possible birth sites of this pulsar, where the path shown intersects the different ^birth values, are highlighted as black, curved sections 
at positions (1) and (2). (b) Greyscale map of the 2D PDF of ^birth an d £kin derived from our MC simulation for PSR J1932+1059 (see 
Eq. I15|) . The two separate regions of high probability in this PDF are the result of this pulsar having crossed the GP twice, an example 
of which is shown in (a), (c) Posterior PDF of ^birth derived from marginalising the 2D PDF over t^ n (see Eq. I16|) . (d) Posterior PDF of 
fkin derived from marginalising the 2D PDF over ^birth ( see Eq. I16|) . The horizontal, dashed line corresponds to n, i.e. the upper limit 
of our simulation; the PDF for larger t^- in is shown for completeness as a grey outline. 



sar's path has intersected with the GP more than 10 times 
(for certain values of v r ). In such cases, although we can 
assign a most probable tkin to these pulsars based on the 
value at the highest peak, we cannot confidently exclude 
ages corresponding to other peaks of roughly equal proba- 
bility to that of the highest. A few of those cases exhibit a 



clear bimodality in the shape of p(tkin), with a peak sepa- 
ration that can be up to 2 orders of magnitude in units of 
t c . Fig. [5] shows a characteristic example of this bimodal- 
ity for the case of PSR J1932+1059, where p(tkin) contains 
two significant peaks near tkin = 0.1r c and 10r c . The na- 
ture of those peaks becomes clearer in the grayscale map 
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of p(vbirth, tkin), where two distinct regions of high proba- 
bility are visible, each corresponding to the pulsar's inter- 
section with the GP (i.e. \z\ ^ 0.1 kpc), separated by the 
interval when the pulsar was travelling outside the plane, 
i.e. \z\ > 0.1 kpc (Fig. [SJd). The intersections are shown in 
Fig. [5ji, at positions (1) and (2), where we have chosen an 
example trajectory corresponding to v r — —250 km s _1 . 
Therefore, when there is an ambiguity in tkin, we report 
only the value corresponding to the highest peak without 
attempting to estimate uncertainties. 

Table [1] shows the resulting tkin values from our MC 
analysis, for all pulsars with convergent solutions within 
tkin ^ Ti, |i7 r | < 500 km s _1 and \z\ ^ 0.1 kpc. All ambigu- 
ous determinations are denoted with an exclamation mark 
next to the value of tkin- 

3.3.2 The reliable sample of kinematic ages 

Despite the difficulties in assigning a single, most-likely kine- 
matic age to several pulsars in our sample, it is still useful 
for all pulsars to report how probable it is that the true 
age — based on kinematics and all our assumptions — lies 
within the intervals demarcated by the characteristic time 
scales t c , T2 and n. These probabilities were calculated for 
each pulsar as 

P3 = p(£kin ^ T c) = / P(£kin)d£kin 
JO 

P2 = P(j c < tkin < T 2 ) = / p(£kin)d£kin (18) 

J T C 

poo 

Pi = P(r 2 < tkin < Tl) = / P(t kin )dt kin 

J T 2 

and are shown in columns 12-14 of Tabled with the highest 
probability of the three shown in boldface. 

In addition, the probabilities p\ , p2 and ps can be used 
to discuss the most-probable range of braking indices for our 
pulsar sample. In Section [5l we investigate the implications 
of the tkin distributions for the distribution of n. 

Our analysis was restricted to tkin ^ n- However, given 
enough time, many of the pulsars in our sample would inter- 
sect with the GP beyond n , effectively providing an infinite 
number of solutions for tkin- For example, the case of PSR 
J 1932+ 1059 (Fig. [5) clearly demonstrates that allowing only 
tkin ^ ti truncates p(tkin) near the position of the most 
probable value of tkin- It is also clear from that plot that 
the PDF exhibits another peak of even higher probability at 
tkin ~ 100t c or 10ti (shown with a grey line). Nevertheless, 
such large departures from the characteristic timescales for 
this pulsar are beyond the scope of the present investigation. 

After excluding all pulsars that did not satisfy the above 
criteria — those marked with an exclamation mark in Ta- 
ble [1] — we were left with 33 pulsars for which tkin was con- 
sidered reliable for the present work. All the derived p(tkin) 
for the selected sample are shown in Fig. [19] Using these pul- 
sars alone, we wish to re-examine the distributions of Fig.[T] 
where the age bins are now defined by the 'most probable 
tkin' instead of r c . Hereafter, when we mention the tkin value 
of a pulsar, we are referring to the most probable value in 
the corresponding p(tkin)- A comparison between the dis- 
tribution of r c , using all 58 PSRs, and that of tkin, using 
the MC values produced by our simulation for the 33 PSRs 



. I . Selected sample (33 PSRs) 

+0.5 

1 .o - <l0 9 (Wy)> = 6 - 2 _o. 6 



Original sample (58 PSRs) 

+0.5 

<log (T c /yr)> = 6.5 _ 1 1 




log (Age) (yr) 

Figure 6. Grey-filled bins: r c PDF of the original sample of 58 
pulsars considered in the herein study of spin- velocity alignment. 
Solid-line steps: tkin PDF generated from all MC values produced 
by our simulation for 33 pulsars with a reliable t^in determination 
(see Table [TJ. Dashed-line steps: r c PDF of 1,500 non-recycled 
pulsars from PSRCAT. The vertical lines of corresponding line 
style indicate the median of each distribution. The medians and 
the 68% CLs are also shown in the inset key. 

that passed our criteria, is shown in Fig. [6] The median 
value for both distributions is also shown in the same figure, 
with <log(t k in/yr)> = 6. 2^06 and <log(r c /yr)> = 6.5l? - .i. 
Although the PDF of tkin appears to be slightly shifted to- 
wards younger pulsar ages compared to the PDF of r c , the 
difference is not significant compared to the la confidence 
intervals of the distributions. 



3.3.3 Very young pulsars 

A notable difference between the distribution of tkin for the 
selected 33 pulsars and that of r c for the original sample 
of 58 pulsars (Fig. [6} is the longer low-age tail of the lat- 
ter, extending below 100 kyr. This is largely due to the 
characteristic ages of the Crab and Vela pulsars. The Crab 
pulsar is missing from the sample of tkin because within 
n = 2r c Crab ln(P C rab/l ms) « 9 kyr, where r c Crab = 1, 240 y, 
it never reaches within 100 pc of the GP. As a result, our 
simulation did not produce any solutions for this pulsar (see 
Fig. Et). The Vela pulsar (r^ ela « 100 kyr) is excluded on 
the same grounds (see Fig. [TJd): it is presently found in the 
plane, but its motion in the z direction is very limited within 
T]^ ela to result in a reliable tkin- 

Amongst the pulsars with a reliable tkin determination 
are PSRs J0538+2817, J1801-2451 and J1952+3252, whose 
ages have been previously estimated from either pulsar tim- 
ing and/or VLBI measurements of their outwards proper 
motion from the respective centres of their parent SNRs 
(Kramer et al. 2003; Ng et al. 2007; Gaensler et al. 2000; 
Zeiger et al. 2008). From that work, the true age of PSR 
J0538+2817 was determined to be 20-40 kyr, roughly an 
order of magnitude younger than its characteristic age. Our 
simulations assign a tkin = 631^473 kyr, for this pulsar, with 
a probability of p(£kin < 40 kyr) « 7%. Furthermore, the 
true age of PSR J1801— 2451 is quite uncertain and could 
be as old as 170 kyr, if it was born in G5.4— 1.2, or as young 
as r c = 15.5 kyr, if the association is false and the charac- 
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Table 1. Age and kinematic properties of 52 pulsars (columns 4-9) — taken from the ATNF catalogue (Manchester et al. 2005) — for 
which we were able to determine a kinematic age, (column 11), in the interval 0— r\ (column 10; see text for definition). The distance 
to these pulsars has been estimated using various methods, e.g. from the pulsar DM and a model of the Galactic electron density, from 
VLBI observations, or from timing parallax. Columns 12-14 show the probability of having t^[ n ^ r c , r c < tkin ^ r 2 and T2 < tkin ^ r i> 
respectively, in the considered age interval. The asymmetric errors on tkin correspond to the 68% CLs around the peak value. In several 
cases, the tki n PDF contains multiple peaks and discontinuities: those cases, denoted with an exclamation mark, result in an ambiguous 
determination of t^in an d the peak value is reported without the corresponding CLs. For those pulsars for which the determination of 
tkin was not ambiguous, the most probable value of the birth period for n = 3, Po,3, and its 68% CLs is provided in the last column (see 
Section [P]) . 
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teristic age is close to the true age. Our simulation yields 
tkin ^ 158 kyr, with a 92% probability that tkin > t c , which 
is consistent with the above range of estimates. Lastly, PSR 
J1952+3252 has been firmly associated with CTB 80 which 
implies a true age of ~ 51 kyr in contrast to its ~ 100- 
kyr characteristic age. Our simulation does not predict such 
a young age, producing instead tkin = 631^233 kyr, with 



p(tkin ^ 96 kyr) = 0. Similarly to the cases of the Crab and 
Vela pulsars, the above comparisons highlight the weakness 
of our kinematic method in predicting the ages of very young 
pulsars, which have not traversed significant lengths through 
the Galaxy, away from their birth locations. In Section ^. 4. 11 
we explore the impact of these young pulsars on the overall 
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Table 1. Continued. 
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(!) Complex PDF with multiple significant peaks and in some cases discontinuities. Only the value corresponding to the PDF maximum 
is tabulated. 



trend of spin-velocity correlation with pulsar age, by con- 
sidering their true ages. 

In conclusion, our selection based on the degree of am- 
biguity in the determination of tkin results in a somewhat 
narrower distribution of ages than that of r c , which is also 
quantitatively revealed by comparing the variances of the 
two distributions: Var(r c ) = 0.9 and Var(tkin) = 0.4. This 
difference is largely due to the two youngest pulsars missing 
from the selected sample: the Crab and Vela pulsars. How- 
ever, the selected sample is representative of kinematic ages 
covering a wide range of pulsar ages, ranging from ~ 0.1 
to 100 Myr; this allows us to use the derived tkin values as 
an independent variable in the investigation of spin-velocity 
alignment as a function of pulsar age, like we did in Section[2] 
with r c . 

3.4 Spin— Velocity Alignment Revisited 

The p(tkin) derived in this work provide a quantitative tool 
for assigning an age to each pulsar, at a certain confidence 
level, that is an alternative to the characteristic age. Hence, 
it is also interesting to re-examine the spin-velocity align- 
ment as a function of pulsar age, as was done in Fig. [2] (see 
Section [2]), in the light of the newly derived kinematic ages 
for the pulsars in our sample. As was mentioned earlier, 
our re- investigation of the spin-velocity alignment included 
only the 33 pulsars of Fig. [19] Similarly to Fig. [T] the tab- 
ulated distributions of ^F, for intervals of tkin instead of r c , 
are shown in Fig. [8] As before, in this plot each cell contains 



the distribution of ^ generated from pulsars in intervals of 
tkin, ranging from 100 to 100,000 kyr; also, the correspond- 
ing KS probabilities of rejecting uniformity are shown with 
each distribution. For comparison, we have overlaid the cor- 
responding distributions of Fig. [T] with dashed lines. Simi- 
larly to Fig. [21 a scatter plot of the degree of spin-velocity 
correlation as a function of tkin is shown in Fig. [9] 

The direct comparison between pKs(tkin) and pks(t c ) 
reveals a very similar behaviour for pulsar ages < 10 Myr, 
with both distributions rejecting that spin and velocity 
are uncorrelated with > 80% confidence. There are indeed 
small differences between the central values of pKs(tkin) and 
Pks(Vc), but these could be partly attributed to the differ- 
ent sample sizes; furthermore, if we take into account the 
68% CLs of pks, both cases are consistent with each other. 
However, above 10 Myr, which roughly equates to tdyn, the 
PDF of ^ is evidently more uniform compared to those for 
younger pulsars; this is also confirmed by the small value 
of PKs(tkin) ~ 0.34. The flattening of the ^ distribution for 
samples containing older pulsars is more clearly displayed in 
the scatter plot of Fig. [9l where a downwards trend for pks 
can already be seen across the 1-10 Myr range — although 
the uncertainties on pks are quite significant. The ^ distri- 
bution for pulsars with 10 < tkin < 100 Myr is evidently 
more uniform than those based on younger pulsars, which is 
confirmed by the low value of pks - However, we should warn 
that the bin corresponding to the oldest kinematic ages con- 
tains only 3 pulsars, which is just below the margin for a reli- 
able report of the KS statistic (Stephens 1970). Interestingly, 



© 0000 RAS, MNRAS 000, 000-000 



Pulsar Spin-Velocity Alignment II 13 



0.15 
0.1 
0.05 


-0.05 
-0.1 

-0.15 
-0.2 

-0.25 



- 


x \v.a%\% - 








(a) 


PSR J0534+2200 " ® 



-1 .5 -1 

y(kpc) 



(b) : 




f' 


















PSRJ0835^510 





-0.6 -0.55 -0.5 -0.45 -0.4 -0.35 -0.3 -0.25 -0.2 

y(kpc) 

Figure 7. Example trajectories of the (a) Crab and (b) Vela 
pulsars, on the y—z plane, between their present positions (circled 
dots) and those in the past, delineated with grey dashed lines, for 
various magnitudes of the radial velocity, v r - Our simulation only 
considered path lengths of up to r\ in the past (black solid lines); 
for the Crab pulsar, t\ is very short and cannot be seen on the 
chosen scale. The investigated range of birth heights in our work 
was l^birthl ^ 100 pc, which is demarcated with black, horizontal, 
dotted lines. 



as can be seen in the same plot, spin-velocity alignment is 
retained at a high confidence (pks(t c > 10 Myr) ~ 98%) for 
pulsars in the same age interval according to r c . Nonetheless, 
it should be noted that of the 12 pulsars with r c ^ 10 Myr 
in the original sample, only PSR J1456— 6843 has a reliable 
tkin m the corresponding age interval; the rest of the pulsars 
have complex p(tkm) an d were therefore rejected. Two pul- 
sars, PSRs J0452-1759 and J1937+2544, have a t kin that 
is roughly an order of magnitude older than their charac- 
teristic ages t c < 10 Myr; these are included in the top age 
interval of the new table. Lastly, since none of the pulsars in 
our sample have tkin ^> 100 Myr, we are not able to exclude 
spin-velocity correlation for very old pulsars. However, given 
the observed trend of pks with tkin, we consider it very likely 
that this is the case. 

In conclusion, our alternative investigation of spin- 
velocity alignment as a function of pulsar age, based on 
pulsar kinematic ages, resulted in an intriguing flattening 
of the ^ distribution beyond 10 Myr. Although only 3 pul- 
sars have tkin > 10 Myr, the overall downward trend in pks 
for tkin < 10 Myr can indicate that the Galactic gravita- 
tional potential plays indeed an important role in smearing 
out the observed spin-velocity alignment for older pulsars. 



3.4-1 Very young pulsars 

There are a number of young pulsars in our sample for which 
we have independent estimates of their true ages, based on 
their SNR associations. Furthermore, the spin orientations 
of most of these young pulsars have been derived from fits to 
pulsar wind nebulae (PWNe) seen in their X-ray images (see 
Table 2 of Ng & Romani 2004). So, given the well- measured 
orientation of these pulsars, it is interesting to discuss their 
impact on the overall picture of Fig. [8] by considering their 
true ages. 

In the previous section, we mentioned the age esti- 
mates for PSRs J0538+2817, J1801-2451 and J1952+3252 
from VLBI measurements of their outwards proper motion 
in their SNRs. In addition, the Crab and Vela pulsars, for 
which we were not able to derive a reliable tkin, have well- 
known true ages: the Crab pulsar has a historical age of 958 
y, whereas the outwards proper motion of the Vela pulsar 
from the centre of its SNR implies an age of 18 =b 9 kyr 
(Aschenbach et al. 1995). Moreover, the age of the Vela pul- 
sar has been independently estimated to be 20-30 kyr, from 
measurements of the pulsar's braking index (n ~ 1.4) and 
assuming a birth period of Po w 20 ms (Lyne et al. 1996). 
Also, the spin directions of these young pulsars from the 
PWNe fits result in small spin- velocity alignment angles, 
with ^ - 5°-12°. The exception is PSR J1801-2451, for 
which the spin-axis orientation has been derived from polar- 
isation measurements and results in ^ = 35° ± 5° (Johnston 
et al. 2007). 

We have calculated the probability of rejecting unifor- 
mity, pks, using only the aforementioned 5 pulsars, using 
the independent estimates of their true age. Apart from the 
age of PSR J1801— 2451, for the rest of the pulsars these 
estimates are significantly younger than both r c and tkin, 
which securely places them in the 0-100 kyr interval. If 
we use the published ^ values for these pulsars, including 
that for PSR J1801— 2451, then uniformity is rejected with 
only pks — 76^34% probability. However, as was stated ear- 
lier, the age of PSR J1801— 2451 is quite uncertain and it 
could be significantly older than 100 Myr (Gaensler & Frail 
2000). On that account, if we exclude this pulsar from the 
above sample, we obtain a much more significant rejection 
of uniformity, at the level of 94^i 6 %. The overall trend of 
spin-velocity correlation with age is shown in Fig. [9l where 
we have also added the corresponding pks value for the 0- 
100 Myr interval from the 4 young pulsars. It can be seen 
from that plot that our general conclusion of a strong spin- 
velocity correlation for 10 Myr that diminishes beyond 
10 Myr is maintained. 



4 DISCUSSION 
4.1 Pulsar Distances 

In this section, we discuss another important factor that can 
potentially affect the above conclusions: the large uncertain- 
ties on pulsar distances. As was mentioned in the Introduc- 
tion, DM-based estimates of pulsar distances are often sig- 
nificantly different from distances derived from SN associa- 
tions, HI absorption or parallax measurements (e.g. Chater- 
jee et al. 2009; Deller et al. 2009; Schnitzeler 2012). Fig. [TUl 
shows a scatter plot of the best available distances for our 



© 0000 RAS, MNRAS 000, 000-000 



14 Noutsos et al 



100 



1,000 



10,000 



100,000 








Pks = 



N =30 



10 20 30 40 

T(deg) 




10 20 30 40 

T(deg) 




10 20 30 40 

T(deg) 



100 




10 20 30 40 

Y (deg) 




10 20 30 40 

Y (deg) 




10 20 30 40 

Y (deg) 



1,000 




Pks = 0-726 




10 20 30 40 

V (deg) 



20 30 

V (deg) 



10,000 




Figure 8. Tabulated distributions of \£ in intervals of tkin, for 33 pulsars of the original sample of Fig. [T] after omitting 19 pulsars 
with poor tkin determination. The details of the statistical analysis are identical with those for Fig. [T] (see corresponding caption). The 
^ distributions shown here can be directly compared with those of Fig. ^ for the corresponding intervals in r c : these are shown with 
dashed, red lines on the same scale. Our analysis did not produce reliable kinematic ages below 100 kyr, so the corresponding interval in 
this plot is empty. A discussion on pulsars with very young ages can be found in Sections 13.31 and 13.41 



sample of pulsars (black bullets), as well as alternative dis- 
tance estimates based on the NE2001 electron-density model 
(grey bullets; Cordes & Lazio 2002); also, where available, 
lower and upper limits on the distance based on HI absorp- 
tion and SNR associations are shown (triangles). It is clear 
from that plot that the differences between distances based 
on the model and those based on alternative methods can 



range from a few percent to several tens of percent; and in 
some cases the difference can be of the order of a few hun- 
dreds percent. It is hence important to examine the impact 
that a wrong estimate on the distance would have on the 
kinematic ages. 

We have reproduced the tkin distributions for all pulsars 
in our sample assuming their distance is 0.5 and 2 times the 
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Figure 9. As in Fig. [2] but for intervals of the kinematic age, 
tkin, an d for the remaining sample of 33 pulsars, after excluding 
pulsars for which tki n was ambiguous. For comparison, we have 
over-plotted the central values of Pks( t c) with light-grey dots, 
assuming the same scale for r c as that of tki n shown. The left- 
most data point, for tkin < 100 kyr, shown with dashed lines, 
corresponds to the distribution of ^ values for the Crab and Vela 
pulsars, and PSRs J0538+2817 and J1952+3252, if these pulsars 
are included in our analysis with their true ages: tcrab — 958 y; 
tvela ~ 25 kyr; tjo538+28i7 ~ 30 kyr; and tji952+3252 ~ 51 kyr 
(see text for references). 



published distance, d; the choice of varying the published 
distance by factor 2 was motivated by the worst-case dis- 
crepancies between the various distance estimates shown in 
Fig. [10] The resulting values of tj^f , for pulsar distances 
equal to 0.5d, and t^m, f° r pulsar distances equal to 2d, are 
shown in Table [2] A visual comparison between all three 
cases, i.e. for pulsar distances equal to d, 0.5d and 2d, is 
provided in the scatter plot of Fig. II Ik . In a few cases, our 
simulation did not converge to a solution within n , when the 
alternative distances were used: the data points correspond- 
ing to those cases are missing from the plot. Furthermore, 
as was done for Table [1] for pulsars with an ambiguous de- 
termination of £km, we report only the kinematic age corre- 
sponding to the maximum of the PDF, and without quoting 
errors; the corresponding data points are also omitted from 
the scatter plot of Fig.[TT] For reference, the same plot shows 
the value of r c , for each pulsar. 

At first glance, Fig. Ilia reveals that for the major- 
ity of pulsars both tjVf an< ^ ^kfn are consistent with the 
values of tki n , within the quoted errors. In order to as- 
sess this apparent similarity more quantitatively, we calcu- 
lated the difference between logt^m an d log^khf > i- e - = 
log tkin - log tkhf and between logtk in and log t£f n , i.e. A 2 d = 
logtidn — l°g^k?n 5 at every MC iteration of our simulation. 
The distributions of A^/ 2 and A 2 d, for all 33 pulsars are 
shown in Fig. II 1 b . It is clear from these distributions that, 
on average, the modification to the pulsar distance does not 
significantly impact on tkin, considering the errors on the 
latter. More specifically, we find that (A d / 2 ) — — 0.01^0 53 ? 
where the quoted 68% interval is comparable to the standard 
deviation of the p(tkin) for the same sample (see Fig. [6}. On 
the other hand, the distribution of A 2c z appears somewhat 
narrower, with (A 2c z) = — O.Ollo'J^ revealing the lesser im- 
pact on tkin by doubling the distance compared to that by 
halving it. 

Overall, we can conclude that the impact of distance 
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Figure 10. The best distance estimates of the 52 pulsars in our 
sample (black), contrasted with those based on the pulsar DM and 
the NE2001 model (gray), and the lower (red) and upper (blue) 
limits based on HI absorption or SNR associations. The index 
number iV on the horizontal axis corresponds to that shown in 
the first column of Table [T] 



errors on the determination of tkin, at least for our pulsar 
sample, can be neglected in comparison to the uncertainties 
arising from the unknown v r and Zbirth- There are 3 excep- 
tions to this in our sample, PSRs J0452-1759, J0538+2817 
and J 1937+2544, for which the difference in tkin is greater 
that la between the different distance assumptions. The rea- 
son for the significant impact of the distance on the kine- 
matic age for these 3 pulsars cannot be simply traced down 
to having special properties: neither their positional (/,&) nor 
their kinematic (vi,Vb) parameters are outliers in our sam- 
ple, which would account for the divide. Clearly, the impact 
of distance errors on estimates of tkin warrants a complex 
description. 



4.2 Pulsar Braking Indices 

4-2.1 In general 

Despite being very difficult to measure, there are 8 pul- 
sars so far with reliable braking indices in the literature 
(Lyne et al. 1993,1996; Middleditch et al. 2006; Livingstone 
et al. 2007; Weltevrede et al. 2011; Espinoza et al. 2011). 
Although the pulsar sample is small, it is clear that none 
of them, in fact, have n — 3, whereas the majority have 
2 < n < 3 with the rest of the values being ~ 1. Further- 
more, if we assume that the PDFs of tkin derived from our 
simulation correspond to those of true pulsar ages, we see 
that for a number of pulsars in our selected sample r c falls 
several a away from the most likely value of tkin- The pos- 
sible reasons for such significant discrepancies can be better 
understood if we examine the dependence of n on tkin and 
Po, as is expressed by Eq. [2] 



n = 1 + 2- 



tkin 



(19) 



As n is non-separable from Po, in this equation, we can only 
solve it numerically for a range of t^i n /r c and Po/P values. 
Fig. IT2l shows n as a function of Po/P, for different values 
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Table 2. Kinematic ages of 33 pulsars with reliable estimates of t^ n , assuming they are at different distances: (Column 4) kinematic 
age assuming the published pulsar distances, d (shown in column 7 of Table [TJ, i.e. logf^ in ; (Column 5) kinematic age assuming half 
the published distance, i.e. logt^hf ; and (Column 6) kinematic age assuming twice the published distance, i.e. \ogt^.f n . Undeterminable 
solutions of tki n , for a particular distance assumption and pulsar are shown as dashes. Pulsars for which the fkin determination is 
ambiguous for the particular distance assumption are accompanied only by the t^in value corresponding to the PDF maximum. Also, for 
reference, Column 3 shows the pulsar characteristic age. A scatter plot of the tabulated data is shown in Fig. [TT] 
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(!) Complex PDF with multiple significant peaks and in some cases discontinuities. Only the value corresponding to the maximum of 
the PDF is tabulated. 
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Figure 11. (a) Top panels: Comparison of the kinematic ages of 
the 33 pulsars in our selected sample (black) with those derived 
using the same analysis but assuming the pulsars are at half their 
published distance (red) and at twice their published distance 
distance (blue). For the sake of comparison, the value of r c for 
each pulsar is also shown as a grey dot with a horizontal bar. The 
index number N on the horizontal axis corresponds to that shown 
in the first column of Table ^ Missing values in this plot imply 
that no solution was found within r\ for the corresponding pulsar, 
using the corresponding distance, (b) Bottom panels: (left) PDF 
of the difference A d / 2 , between the log tkin values corresponding 
to distance d and to d/2, for the selected 33 pulsars, and (right) 
of the difference A 2 d between the log tki n values corresponding to 
distance d and to 2d. The value of the median of A d / 2 an d 
(black vertical lines) and their 68% CLs (dashed vertical lines) 
are shown in the inset keys. 




Figure 12. Pulsar braking index, n, as a function of birth period, 
Po (normalised by the present period, P), for a range of tkin Ac 
values. 



In contrast, pulsars with £ kin ^> r c , like PSR 
J0922+0638, pose an interesting problem, as the discrep- 
ancy between tkin and r c cannot be reconciled under n = 3, 
for any value of Po/P (N.B. pure magnetic dipole braking 
only allows tkin ^ t c ). For those cases, we are forced to ac- 
cept that n < 3, if the pulsar was born close to the GP. In 
our selected sample, there are 7 pulsars for which we can 
confidently exclude that tkin < t c under the assumptions of 
our simulation. Critically, as we saw in Section T4. II an error 
on the assumed pulsar distance could modify p(tkin) to a 
certain degree. However, this was found to be a mitigating 
factor only in one case, that of PSR J0738— 4042, whose tkin 
could be made consistent with r c if the pulsar is assumed to 
be 50% closer to the Sun (see Table [2} . For the rest of the 
pulsars, a modification to their distance by up to a factor of 
2 cannot reconcile the fact that p(tkin) is inconsistent with 
pure magnetic dipole braking. 



of tkin/^c. The most evident feature of n(Po/P) in this plot 
is that n approaches its maximum value of 1 + 2(r c /tki n ) as 
Po/P -> 0, for any t k i n /T c . 

Moreover, as was also discussed in Section [3T] it is clear 
that for a constant n, a larger pulsar age implies a smaller 
Po and vice versa: i.e. a pulsar born slow-spinning would 
spin down to its present spin period, under a constant n, 
faster than the same pulsar born fast-spinning. Hence, pul- 
sars like PSR J1955+5059, for which tkin ~ 0.1r c , may have 
been born with Po ~ 482 ms and spun down with a con- 
stant n — 3 or possibly with Po <C P and spun down with 
n = 1 + 2r c /tkin ~ 17, or perhaps, what is a more likely 
scenario, with another combination of n and Po that lies 
between those values of n and Po/P, for t^i n /r c ~ 0.1. Al- 
though, it should be noted that the value of n for such small 
age fractions remains roughly constant and equal to its max- 
imum allowed value for up to Po/P ~ 0.95, as can already 
be inferred from Fig. 1121 Hence, unless we are willing to ac- 
cept extreme values of n for PSR J1955+5059, the results 
of our simulation suggest that it was very likely born with 
Po - P. 



4-2.2 Application to the selected sample 

The parametric curves of n(Po/P) in Fig. [12] give us a gen- 
eral idea of the expected distribution of n for pulsars of 
different ages, born spinning slowly or spinning fast. In par- 
ticular, for Po/P ~ Owe have n ~ l+2(r c /tkin) > 1, with the 
upper limit restricted only by the minimum value of t^i n /r c 
— e.g. the lower la CL on t kin /r c for PSR J0837-4135, be- 
ing « 0.02, would yield n ~ 100. On the other hand, for 
Po/P 7^ 0, the range of values for n is unrestricted, and can 
even be negative: the lower bound on n only comes from 
the largest Po/P we are willing to accept. In those cases, 
small values of t^i n /r c (e.g. < 1), as in the case of PSR 
J1955+5059, would yield n < 1 + 2r c /t kin only for P ~ P, 
whereas pulsars with large tki n /T c , like PSR J0922+0638, 
would yield n < for all but the smallest values of Po/P- 

One of the interesting open questions regarding pulsar 
evolution is the dominant spin-down mechanism responsible 
for the bulk movement of pulsars on the P—P diagram. The 
small number of measured braking indices does not allow 
for a statistical analysis. However, our simulation allows us 
to calculate the distribution of n from p(tkin), for each of 
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Figure 13. PDFs of the braking index, n, calculated from the 
PDFs of t kin and Eq. [19] for the 33 selected pulsars of Table [2] 
The distributions shown correspond to (from top to bottom) 
Po/P = 0,0.01,0.1 and 0.5. The vertical dotted lines demarcate 
the intervals between n = 0, 1, 2 and 3, where we also show the 
percentage of values that fall into each interval. 



the selected 33 pulsars, using Eq. [19] and assuming a value 
for Po/P. The aggregate distribution of n for all 33 pulsars 
assuming P /P = 0, 0.01, 0.1 and 0.5 is shown if Fig. ED As 
expected, assuming that all our pulsars were born spinning 
infinitely fast (Po/P = 0) results in all values of n being 
distributed above 1, with a noticeable peak near n = 1. The 
reason for the accumulation of values in the bottom bin of 
that distribution is that all solutions with £ki n ^> T2 will be 
binned near n = 1, since slower spin-down is not allowed. In 
other words, the distribution of n for Po/P = is populated 
with values that are in the majority upper limits to the true 
values of the braking index. If we allow finite values of Po, it 
can be seen in the rest of the distributions that a significant 
fraction of solutions previously corresponding to n « 1, is 
now distributed between n = and n = 1. Finally, in the 
more extreme case of Po/P = 0.5, the PDF of n is practically 
uniform across a wide range of n and even allows for n < 0. 

In conclusion, based on the distributions of n for our 
sample, under different assumptions for Po/P, it is difficult 
to determine the most likely range of braking indices, inde- 
pendently of Po - Indeed, if we assume that all 33 pulsars were 
born with Po/P < 0.01, then n ~ 1-2 seems the most likely. 
However, we cannot exclude the possibility of pulsars being 
born with spin periods much nearer their measured value; 
in which case, as can be seen in Fig. [13] spin-down models 
with < n < 3 appear to be an equally likely scenario. 




100 



200 

P 03 (ms) 



300 



Figure 14. The birth-period PDF, p(Po,3), calculated from the 
kinematic-age distribution of PSR J0454+5543, p(t^ n ), assuming 
n = 3 (pure magnetic-dipole braking). The solid line corresponds 
to the peak of the PDF and the dashed and dotted lines, to the 
68% and 95% confidence limits, respectively. The rightmost, dash- 
dotted line marks the presently measured period for this pulsar. 



4.3 Pulsar Birth Periods 

4-3.1 In general 

In addition to the determination of n from tkin, assuming 
different values for Po/P, one can invert the problem and 
determine Po for different spin-down models, assuming a 
value for the braking index, n: i.e. by setting ttme = tkin 
and solving Eq. [2] for Po, 



l)tkh 



2r c 



(20) 



where P is the observed pulsar period. Note that the maxi- 
mum allowed value of tkin for a given n, for which Po remains 
a real number in this equation is r n = 2r c /R(n — 1); here, 
R(n — 1) is the ramp function, which is equal to n — 1 for 
n ^ 1 and everywhere else. This is the generalised form 
of the characteristic age for different n, which assumes that 
Po/P = 0. The reader will notice that, according to the 
above definition of r n , the characteristic age for n = 1 is 
infinity, which is true if Po = is assumed; but n is finite as 
we have defined it in Eq.0 The allowed range of t^i n /r c for 
different spin-down models is schematically shown in Fig. 1151 
Assuming a pure magnetic-dipole spin-down, i.e. n = 3, 
Eq. 1201 becomes 



Po,3=P\ll 



^kin 



(21) 



In this section, when appropriate, we use two subscript in- 
dices for P: the first index denotes the time at which the 
period is calculated (e.g. at birth, t = 0), and the sec- 
ond element, the assumed braking index in that calculation 
(e.g. n = 3). 

We have used Eq. [21] to calculate the distributions of 
Po,3 for each of the 27 from the 33 selected pulsars, for which 
P3 = p(tkin ^ t c ) > 0. Fig. 1 141 shows the distribution of Po,3 
for PSR J0454+5543. The value of P corresponding to the 
maximum of the PDF for each of the 27 pulsars is shown in 
column 15 of Table [T] 

In several cases, the calculation of the PDF of Po,3, from 
Eq. 1211 excludes a significant fraction of p(tkin)- An exam- 
ple that highlights this fact is that of PSR J1453+6413, for 
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Figure 15. Pulsar spin period at birth, Po (normalised by the 
present period, P) as a function of the kinematic age, tkin (nor- 
malised by r c ), for a range of braking index values. 



which Po,3 can be calculated from only 15% of the values 
produced in our simulation. Similarly, for n — 2, Po,2 can 
be calculated for about 74% of those values; and, finally, 
Po,2 can be calculated for all values of tkin- The shape of 
these PDFs depends directly on n (Eq. EDI), but 

since dif- 
ferent choices of n also allow different ranges of tkin, n also 
influences the PDFs indirectly: 

(a) The direct dependence on n can be estimated by calcu- 
lating the first derivative of Eq. [20] with respect to n. The 
fractional change of Po, SPo/Po, as a function of the change 
in n, Sn, is then 



SPp 
Po 



2(1 



tkin 



Sn 



Sn 



if n : 



(22) 

It can be shown that SPq/Pq always has the opposite sign 
to Sn. Therefore, assuming a smaller value of n results in 
larger values of Po for the same tkin- 

(b) Further to the above, if we assume a small change of 
the braking index by Sn <C n, then the range of tkin that 
contributes to p(Po) will change by 



Sr n 



R(n + 6n-l) R(n - 1) 



(23) 



So, an increase in n shrinks the allowed range of tkin and 
a decrease, expands it. The effect of Sn on the distribution 
of Po depends on the value of tkin relative to r m , where 
m — max (n, n + Sn). In general, the distribution of Po will 



be modified according to the following: 

p(^ T n ) 



p[Po,n+<5n(tkin) 



T n +Sn) 



p[Po,n(t k in)] if tkin ^ T m 



H(T n+S n ~ tkinMPo 

if tkin > Try 



where p(^ r n ) = / p (tkin) catkin- Although it is possible, as 
Jo 

was shown earlier, to express the change in P for tkin ^ T m 
analytically, the form of p(Po) for tkin > T m depends entirely 
on p(tkin) and is therefore unique to every pulsar. 

Bearing the above in mind, in Fig. [16] we show two ex- 
amples that highlight the direct and indirect effect of the 
choice of n on the distribution of Po. In those examples, we 
show the distributions p(Po,i), p(Po,2) and p(Po,3) for PSRs 
J1735— 0724 and J1453— 6413, as they were derived from our 
analysis assuming 3 different values for n. Firstly, in Fig.fTSa 
it can be seen that the PDF of Po is progressively shifted to 
larger values of Po with decreasing n; this is also the case for 
the most probable values of Po in those distributions. This 
can be understood by looking at the p(tkin) of this pulsar 
in Fig. 1191 which lies entirely below tkin = t c ; hence, in this 
case, the effect of Sn is solely due to the direct change of 
Po as is predicted by Eq. [22] On the other hand, Fig. [T6b 
shows that for PSR J1453— 6413 the different choices of n 
yield distributions with very different shapes. Also, the most 
probable values of Po in these distributions do not change 
monotonically with n. This can be explained by the fact that 
the 3 different choices of n sample significantly different frac- 
tions of the p(tkin) of this pulsar (Fig. [T9]) . Therefore, the 
shape of the PDFs of Po is determined both by the direct 
change of Po as a function of n and by the indirect influence 
of n, as is described by Eq. [24] 

4-3.2 Application to the selected sample 

In the previous paragraphs, we highlighted the fact that the 
shape of p(Po, n ) is sensitive to the choice of n. Keeping this 
in mind, we examined the aggregate distribution of the birth 
period by binning all Po,3 values of the 27 pulsars, as we did 
for n in the previous section. The normalised PDF, p(Po^), 
is shown in Fig. 117a . 

The derived distribution of the birth periods of 27 non- 
recycled pulsars from our work can be compared with the 
corresponding distribution that was published recently by 
Popov & Turolla (2012): in that work, the authors inferred 
the birth periods of 30 pulsars based on their associated- 
SNR ages. As in our case, their work assumed standard 
magnetic-dipole braking (n = 3) and concluded that a Gaus- 
sian distribution with a mean and standard deviation of 100 
ms — truncated at Po = — is a good description of the 
derived Po values. Their distribution is shown with a dotted 
line in Fig. 117a . The most apparent difference between the 
two distributions is the much broader range of Po predicted 
by our simulation as opposed to the narrow distribution of 
Popov & Turolla. This dissimilarity is confirmed by the two- 
sample KS test, which yields ~ probability of these distri- 
butions coming from the same underlying distribution. Also, 
apart from the main peak at Po ^ 60 ms, the derived distri- 
bution from our simulations exhibits a second long-period 



© 0000 RAS, MNRAS 000, 000-000 



20 Noutsos et al 




P (ms) 

Figure 16. Birth-period PDFs for (a) PSR J1735-0724 and 
(b) PSR J1453— 6413, for different values of the braking index, 
i.e. n = 1, 2 and 3, as indicated. 



component for Po > 600 ms. The latter is mainly due to 
the calculated values for PSRs J0630-2834, J0837-4135, 
J0837+0610, J1136+1551, J1740+1311 and J2157+4017, 
which have P > 700 ms and for which the probability of 
having tkin <C r c is high (see Table [T]) . Critically, these 6 
pulsars are old and without SNR associations, and so they 
were not included in the sample of Popov and Turolla. 

A number of SN-kick models predict a connection be- 
tween Po and ^ (see e.g. Ng & Romani 2007). Accord- 
ing to these population-synthesis models, pulsars born with 
Po < 20 ms and those born with Po > 100 ms tend towards 
alignment. For the short-period pulsars, the fast spin of the 
proto-neutron star tends to average out the kick directions, 
with the net kick being along the spin axis. Conversely, for 
the long-period pulsars, the proto-neutron star must have 
been spinning slowly at the start and would have been kicked 
radially along its spin axis, as to not have been spun up and 
acquire an orthogonal Vbirth- According to the above simu- 
lations, pulsars with intermediate birth periods of 40-60 ms 
are expected to have ^ ^> 0°: this is because these pulsars 
were spinning slow enough that the bulk of the initial kick 
was delivered at large angles, almost perpendicularly to the 
spin direction. 

Our simulation gives us the opportunity to test the 
above predictions, using the Po,3 distributions of the 27 
pulsars. We followed the same procedure that was used in 
Figs. [2] and [9] wherein we calculated the corresponding KS 
probabilities of rejecting the uniform ^ distribution for a 
number of Po bins. In order to roughly match the scale over 
which changes in spin-velocity alignment are predicted ac- 
cording to the SN-kick models, we binned the values of Po 



into intervals of 50 ms, for Po < 500 ms, and into inter- 
vals of 250 ms, for the rest of the Po range: the upper limit 
for the finer binning was chosen to lie roughly on the bor- 
der between the short-period and long-period components 
of p(Po,3) fFig.[T7a). The la confidence intervals of the Po,3 
distributions for each pulsar are significant, i.e. comparable 
to the size of the selected bins. And so, instead of assigning 
the most probable value of Po to each pulsar from its corre- 
sponding distribution, we randomly selected a large number 
of Po values from each pulsar's distribution: for each ran- 
domly selected set of 27 Po values, we calculated pks only 
for those pulsars whose Po value fell into the investigated 
Po bin. The final value of pks assigned to each bin corre- 
sponded to the most-probable value in the distribution of 
all Po values in the bin considered. Also, the 68% CLs were 
calculated around the most-probable value. In this way, we 
incorporated p(Po,3) into the final value of pks- Fig. 117b 
shows a scatter plot of the values of pks across the binned 
Po range. It is clear from that plot that the small number of 
available pulsars results in only a few pulsars being consid- 
ered per MC iteration, in each bin. Consequently, the 68% 
confidence intervals dominate over any possible dependence 
of pks on P . 

In summary, our efforts to extract information on pul- 
sar birth periods based on kinematic simulations and an 
assumption of their spin-down were generally plagued by 
low statistics. Any tests against model predictions demand 
a relatively high resolution in Po space (i.e. SPo ~ 20 ms), 
which our limited sample could not fulfil — allowing in the 
vast majority of cases for only 4 pulsars to be considered 
per 50-ms-wide Po bin. As in the case of pulsar ages, only a 
larger sample of suitable pulsars would help overcome these 
difficulties. 



5 SUMMARY & CONCLUSIONS 

In this paper we investigated whether the high degree of 
correlation between the spin and velocity axes of pulsars 
depends on pulsar age. In the first part of our study, we 
used the characteristic ages (r c ) of 58 non-recycled pulsars 
and showed that the strong correlation does not vanish for 
pulsars with characteristic ages above ^10 Myr, while we 
would expect any nascent correlation to have been washed 
out due to movement through the Galactic potential. This 
expectation is justified by the distribution of time intervals, 
At, during which the velocity direction of pulsars is signifi- 
cantly modified due to the Galactic gravitational potential. 
The distribution of At was calculated for the pulsars in this 
paper, for a change in their velocity direction by AO = 25° 
due to the galactic potential of Paczynski (see Fig. I18[) . The 
time interval corresponding to the peak of the distribution 
is At = 18^5 8 Myr, which is close to an estimate of the 
dynamical time for the Milky Way (see Introduction). 

The conclusion of the first part of this work was based 
on pulsar characteristic ages, which are often unreliable: 
these are derived assuming that all pulsars are born spin- 
ning infinitely fast and that they spin down under a con- 
stant braking index of 3, by emitting pure magnetic-dipole 
radiation. To reach a more conclusive statement about the 
relationship between spin-velocity alignment and pulsar age, 
in the next part we considered an alternative, more accurate 
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Figure 17. (a) Birth-period distribution, p(Po), derived from a 
large number of simulated values of tkin and Eq. 1211 In this plot, 
the solid line corresponds to the most-probable value of Po, and 
the dashed lines, to the 68% confidence limits around that value, 
respectively. In addition, the distribution of 30 birth periods of 
pulsars associated with SNRs, from the recent work of Popov & 
Turolla (2012) is shown with a dotted line, (b) Dependence of the 
probability of rejecting the uniform ^ distribution, under the KS 
test, on the pulsar birth period, Po, calculated for n = 3. The 
plot details are given in Section 14.31 Both plots are based on the 
selected sample of 27 pulsars, after excluding those with ambigu- 
ous kinematic-age determinations and for which the calculation 
of Pq was possible under the assumption of n = 3. 
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Figure 18. Distribution of the minimum time interval, At, be- 
tween the current position of the pulsars of this paper and their 
past position when their 3D velocity direction differed by 25° 
from its present direction: i.e. A6 = 25°. The value of At corre- 
sponding to the peak of the distribution (solid line) and the 68% 
CLs (dashed lines) are shown in the top right corner of the figure. 



measure of a pulsar's age, based on the observed position, 
distance and proper motion of the pulsar: the kinematic age, 
tkin- The latter is defined as the travel time between the lo- 
cation of the pulsar at birth and that at present, through 
the Galactic gravitational potential. We have calculated the 
most-likely values of tkin for 52 of the pulsars in our sample, 



under the reasonable assumption that they were born some- 
where within 100 pc of the Galactic mid-plane (e.g. Cordes & 
Chernoff 1998; Arzoumanian, Chernoff & Cordes 2002) and 
that they have travelled to the currently observed position 
through the gravitational potential of Paczynski (1990). The 
unknown radial velocities and birth heights were assumed 
to follow the distributions of Hobbs et al. (2005) and Reed 
(2000), respectively, which have been derived from a large 
sample of pulsar proper motions and photographic magni- 
tudes of OB stars. For each pulsar, our simulation was lim- 
ited to a maximum age that we denned as the amount of time 
required for the pulsar to spin down to its present period, 
assuming it was born with a period of 1 ms and spins down 
with a braking index of 1: for the vast majority of pulsars 
in our sample, this limit was at least an order of magnitude 
older than their characteristic age. Our simulation resulted 
in 52 pulsars producing at least one intersection with the 
Galactic plane within the explored time interval. This is the 
largest sample of pulsar ages estimated based on kinematics 
to date. 

Following the £ki n estimates, we tested the degree of 
spin- velocity alignment again, this time as a function of the 
more accurate kinematic age. To eliminate systematic un- 
certainties, caused by ambiguous age determinations due to 
multiple intersections with the Galactic plane, we selected 
only the 33 most reliable tki n : these were selected based on 
the complexity of the derived probability density functions 
of tkin- Unlike the case when r c was used as the pulsar age, 
sorting the selected pulsars based on intervals in tkin re- 
sulted in a diminishing degree of spin-velocity correlation 
with increasing age. Pulsars with ages greater than ~ 10 
Myr, whose 3D velocities have almost certainly been signifi- 
cantly altered by the gravitational pull of the Galaxy, yielded 
no correlation. Hence, the present work — together with our 
preceding study of spin-velocity alignment — suggests that 
the orientations of the spin and velocity axes in young pul- 
sars are correlated, while this is no longer observable for 
older pulsars due to the effect of the Galactic gravitational 
potential. 

Our results are tantalising, but must be strengthened 
with more data: e.g. the interval corresponding to the old- 
est kinematic ages in our sample contained only 3 pulsars. 
Proper motion measurements of old, non-recycled pulsars 
are necessary to improve statistics in that interval. In ad- 
dition, the spin-axis directions of these pulsars can be de- 
termined through polarisation measurements. Currently, a 
number of promising VLBI campaigns, which will result in 
hundreds of new pulsar proper motions, are in operation or 
being proposed (PSR7T, Deller et al. 2011; ell, Vlemmings et 
al., eMerlin Legacy Programme). For many of these pulsars, 
these campaigns will also provide accurate distance mea- 
surements through parallax measurements, independent of 
electron-density modelling. Our work has shown that, for the 
vast majority of pulsars in our sample, a change in the dis- 
tance by as much as a factor 2 (typical worst-case scenario 
for NE2001) does not affect the kinematic-age estimates, 
which remain consistent within la. However, for many pul- 
sars the error on their distance derived from their DM can 
be larger still. Our simulations have shown that the effect 
of distance errors on the kinematic age in our simulations is 
complex — and varies significantly with pulsar position and 
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velocity — so the availability of accurate distances will help 
eliminate such systematic uncertainties. 

Upon the completion of the above VLBI campaigns, 
we will be able to combine new and existing measurements 
and improve the statistics. This will help answer many open 
questions about the nature of spin-velocity alignment, like 
the preference of SN kicks towards alignment or orthogonal- 
ity; it will also allow for a distinction between the degree 
of alignment for pulsars born spinning fast and those born 
spinning slow. Finally, with a larger data set we can explore 
finer age intervals, which will provide a better handle on the 
characteristic time scale of the observability of spin-velocity 
alignment. 

Our work has shown that it is possible to systematically 
estimate pulsar ages independently of the choice of a spin- 
down model, while taking into account the kinematic proper- 
ties of individual pulsars. In this work, we used this method 
to study the evolution of pulsar spin-velocity alignment with 
age, which required a pulsar sample with available absolute- 
polarisation information. Nevertheless, this method can be 
used independently with the hundreds of available pulsar 
proper motions to re-evaluate the ages of pulsars; this is in- 
valuable information that can be used, in future, to study, 
amongst other things, the evolution of pulsars on the P-P 
diagram (Fauchere-Giguere & Kaspi 2006; Keane & Kramer 
2008). In this paper, we performed a preliminary statisti- 
cal interpretation of the distributions of the birth periods 
and braking indices for our pulsars, based on the ratios be- 
tween their characteristic and kinematic ages. In particular, 
we translated the distributions of kinematic age into dis- 
tributions of braking index, under different assumptions for 
the value of the birth period. Under the usual assumption 
of Po/P = 0, we find that there is a clear preference to- 
wards n = 1: however, this is a biased result caused by 
the restriction to n ^ 1 under this assumption. Allowing 
our pulsars to be born spinning slower results in distribu- 
tions of n that are noticeably more uniform, thus eliminat- 
ing any preference towards a particular spin-down model. 
Inversely, assuming pure dipole magnetic braking, we calcu- 
lated the birth-period distributions of 27 of the selected 33 
pulsars from their tkin distributions. The aggregate distri- 
bution of Po from this work is significantly wider compared 
to the birth-period distribution of SNR-associated pulsars. 
This can be explained by the very long birth periods (Po ~ 1 
s; comparable to presently observed periods) of several pul- 
sars for which our simulation yields tkin <C r c . 

In the near future, a large number pulsars will be dis- 
covered in all-sky pulsar surveys, like those currently un- 
derway with the Low Frequency Array (LOFAR; Stappers 
et al. 2011), and the Parkes and Effelsberg radio-telescopes 
(Keith et al. 2010; Ng 2011). Many of these pulsars will be 
discovered at low observing frequencies (~ 100 MHz) and 
will thus be nearby because of signal-scattering limitations 
imposed by the interstellar medium over long distances (van 
Leeuwen & Stappers 2010). Hence, it will be possible to ob- 
tain parallax measurements, as was mentioned earlier, and 
estimate their proper motions and distances. Also, many 
more pulsars will be discovered at high latitudes, where scat- 
tering is less detrimental. The kinematic ages of those pul- 
sars will be more reliably determined by methods such as 
that presented in our work. It will then be possible to con- 
duct a study of the distributions of pulsar braking indices 



using large-number statistics, which was not available at the 
time of this work. Crucially, follow-up polarisation measure- 
ments of all the discovered pulsars will yield information 
on their spin-axis directions. The aggregate sample of hun- 
dreds of spin and velocity directions of known pulsars and 
those discovered in the above surveys will provide an un- 
precedented handle on the intriguing phenomenon of pulsar 
spin-velocity alignment. 

Ultimately, the planned Square Kilometer Array (SKA) 
will discover tens of thousands of Galactic pulsars and will 
have the capability to directly measure proper motions for 
several thousands of them. This will help us address many 
open questions regarding the kinematic properties of pul- 
sars, like which is the birth- velocity distribution of pulsars 
and what is the probability of spin-velocity alignment dur- 
ing a core-collapse supernova kick, for different pulsar prop- 
erties at birth. 
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